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ABSTRACT 

Data selection and the determination of systematic uncertainties in the spectroscopic measurements 
of neutron star radii from thermonuclear X-ray bursts have been the subject of numerous recent 
studies. In one approach, the uncertainties and outliers were determined by a data-driven Bayesian 
mixture model, whereas in a second approach, data selection was performed by requiring that the 
observations follow theoretical expectations. We show here that, due to inherent limitations in the 
data, the theoretically expected trends are not discernible in the majority of X-ray bursts even if they 
are present. Therefore, the proposed theoretical selection criteria are not practical with the current 
data for distinguishing clean data sets from outliers. Furthermore, when the data limitations are 
not taken into account, the theoretically motivated approach selects a small subset of bursts with 
properties that are in fact inconsistent with the underlying assumptions of the method. We conclude 
that the data-driven selection methods do not suffer from the limitations of this theoretically motivated 
one. 

Subject headings: dense matter — equation of state — starsmeutron — X-rays:stars — X-rays:bursts 
— X-rays:binaries 


1. INTRODUCTION 

Thermonuclear flashes on neutron stars have been 
used in the past decade to perform spectroscopic 
measurements of neutron star radii and masses (e.g., 
Maiczyna et al.l [ 200,^ lOzell l2006l lOzel et ^ l2009l 
2011)1: IPoutanen et all 12014) . The rich burst dataset 
(|Gallowav et al.ll2008D from the Rossi X-ray Timing Ex¬ 
plorer (RXTE) has not only allowed the measurement 
of the macroscopic properties of half a dozen of neutron 
stars but also initiated detailed studies of systematic un¬ 
certainties in these measurements. Even though the ma¬ 
jority of sour ces showed bursting b ehavior that is highly 
reproducible (iGiiver et al.ll2012allb[ ) , a subset shows burst 
properties and evolution that are more complex. 

A number of studies explored different approaches 
to identifying outliers in the burst samples of differ¬ 
ent sources that conta minate the statistica l inferences 
in the measurements. In lGiiver et al.l (l2012al lbl). we used 
a data-driven approach that employs a Gaussian mix¬ 
ture Bayesian inference to reject outliers. In particular, 
we looked at the cooling tails of X-ray bursts in the flux- 
temperature diagram, which would have yielded identical 
tracks among bursts of the same source in the absence 
of any astrophysical complexities. The degree of scatter 
we observed in the cooling tracks allowed us to measure 
the level of systematic uncertainty in the measurements, 
e.g., due to obscuration, reflection off the accretion disk, 
or uneven burning on the stellar surface, as well as to 
remove a small percentage of cooling tracks that clearly 
did not behave like the majority. 

In an alternate appro ach, Poutanen et al. (2014; see 
also iKaiava et al.l 12014) developed a theoretically mo¬ 
tivated method to select bursts. They used the burst- 
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iriK ne utron star atmosphere models of ISuleimanov et alJ 
(|2012l) to calculate the evolution of the blackbody nor¬ 
malization and chose only those bursts that followed 
their predicted trends. When applied to the source 
4U 1608—52, they found that this approach selects a very 
small fraction of bursts as acceptable ones and leads to 
substantially larger inferred neutron star radii. 

In this paper, we show t hat the burs t s select ed by 
the data-driven approach of iGiiver et al.l (l2012allbl ) are 
not in conflict wi t h the theoretical expectations of the 
ISuleimanov et all (l2012| ) atmosphere models. The theo¬ 
retical models predict that the blackbody normalization 
evolves rapidly in a narrow range of fluxes close to the 
Eddington limit. However, in most bursts, the instru¬ 
ment limitations of RXTE do not allow resolving this 
fast predicted evolution. When this limitation is not 
ta ken into accoun t in bu rst selection, as was the case 
in IPoutanen et al.l (|2014[) , this procedure biases the se¬ 
lection toward a peculiar set of bursts. We further show 
that this subset of bursts are not photospheric radius ex¬ 
pansion bursts, which is a re quirement in the theore tical 
motivated sel e ction criteria of lPoutanen et al.l (12014) and 
IKaiava et al.l (|2014D . Because of both its practical lim¬ 
itations and its biased results, we argue that this data 
selection method does not lead to reliable radius mea¬ 
surements. 

2. APPLYING THE THEORETICAL SELECTION CRITERIA 
TO RXTE BURSTS 

In theoretical models of bursting neutron star atmo¬ 
spheres, the deviation of the spectrum from a blackbody 
is typically quantified in terms of the color correction 
factor 

/c = Tc/Teff, (1) 

where Tc is the color temperature obtained from the spec¬ 
tral fits and is the effective temperature of the atmo¬ 
sphere. The color correction factor depends on the effec¬ 
tive temperature of the atmosphere, its composition, and 
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of the burst luminosity divided by its maximum value. The verti¬ 
cal grey band corresponds to a 10% flux change, which cannot be 
resolved for the vast majority of radius expansion bursts observed 
with RXTE (see Figure [Sj. Even though the models predict color 
correction factors as large as 1.9, the time resolution used in the 
analysis restricts the maximum effective color correction factor to 
< 1 - 6 - 
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Fig. 2.— The evolution of the X-ray flux (blue) and of the 
apparent angular size (red) during burst #21 from EXO 1745—52. 
In this case, the X-ray flux evolves by 18% between the two time 
bins adjacent to the touchdown point. 

the effective gravitational acceleration, i.e., how close the 
radiation flux is to the Eddington critical value. Because 
/c 1 in all relevant conditions, the normalization of the 
blackbody, defined as 

K = {Rbb/D)^ = ^, ( 2 ) 

where F is the flux of the thermal emission and as is the 
Stefan-Boltzmann constant, is not equal to the apparent 
angular size of the neutron star, {Rapp/D)^, but rather 
to 

K = f-\Rapp/Df. (3) 

In Figure [U we plot a typical evolution of the color 
correction factor for tw o of the models described in 
ISuleimanov et all (l2012h . It is evident in this figure that 
most of the evolution of the color correction factor occurs 
within 10% of the maximum flux (i.e., the Eddington flux 
corrected for temperature effects). 

For a given neutron star, {R^pp/D)^ is fixed, and there¬ 
fore, the blackbody normalization is expected to scale 
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Fig. 3. — (Top) The cumulative distribution of the fractional 
burst flux evolution between the two 0.25 s time bins that are ad¬ 
jacent to each touchdown point, for all the photospheric radius ex¬ 
pansion bursts used in lOzel et al.l l|2Q15l) for radius measurements. 
(Bottom) The cumulative distribution of the fractional change in 
the blackbody normalization during the 0.25 s time bin prior to 
each touchdown point, for all the photospheric radius expansion 
bursts used in the lOzel et al.l ll2015f) study. In the majority of 
bursts, the flux evolves by 5-20% and the blackbody normalization 
by 20-80% around the touchdown point, in a way that cannot be 
resolved with current observations. 

as f^. Ideally, this expected evolution of the black- 
body normalization can be used as a selec tion criterium 
(iPoutanen et al.l I20l4 iKaiava et al.l l2014ll as well as a 
way of measuring a combination of the stellar mass an d 
radius (jMaiczvna fc Madeill2005t IPoutanen et al.ll2014ri . 

Unfortunately, a number of complexities make the ap¬ 
plication of this selection procedure to RXTE data im¬ 
practical. First, this method requires the selecting bursts 
that reach or exceed the Eddington limit, which is usu¬ 
ally achieved by searching for evidence for photosperic 
radius expansion (PRE) episodes. Second, it requires 
accurately pinpointing the time of touchdown at the end 
of the PRE event and using only the data after this point 
to measure the evolution of the blackbody normalization. 
Finally, it requires obtaining a large number of counts at 
short integration times in order to resolve the evolution 
of the blackbody normalization with the rapidly decreas¬ 
ing flux. 

We ha ve discussed the crite ria to select bona fide PRE 
bursts in iGuver et al.l (|2012all and will revisit this point 
in the next section. Because of the burst countrates and 
instrumental limitations, burst spectra with RXTE are 
extracted using ^ 0.25 s integrations. Using a typical 
PRE burst from EXO 1745—248, we show in Figure [2] 
how this time binning affects both the time localization of 
the touchdown point and our ability to resolve the rapid 
flux evolution after it. In this example, within one 0.25 s 
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Fig. 4.— The inferred evolution of the color correction factor for 
all the bursts from 4U 1608—52 (modulo an arbitrary normaliza¬ 
tion), as a function of (top) the measured X-ray flux and (bottom) 
the color temperature of the spectrum. Different colors correspond 
to different bursts. In the bottom panel, the evolution of the color 
correction factor from a He-atmosphere model of Suleimanov et 
al. (2012) is overplotted for comparison. Plotting the blackbody 
normalization as a function of the color temperature reduces the 
scatter and allows the theoretically expected correlation to start 
becoming apparent. 

time bin, the flux to the right of the touchdown point 
evolves by 17%, while the photospheric radius to the left 
of the touchdown point evolves by 63%. As can be seen 
from eq. ([3]) and taking into account the dependence of 
/c on the flux (see FigH]), these two changes have op¬ 
posite effects on the blackbody normalization, effectively 
canceling each other out. 

This behavior is not unique to the example shown in 
Figure [2j Indeed, an analysis of all the PRE bursts 
for the six sourc es used in the radius measurements of 
lOzel et all (|2015f) shows a similar level of evolution. We 
show this in Figure [3l where we plot in the top panel the 
cumulative distribution of the fractional change in the 
flux at the two 0.25 s time bins that are adjacent to the 
touchdown point. In the bottom panel, we also show the 
fractional change in the blackbody normalization during 
the 0.25 s time bin prior to the touchdown point. The 
typical values of the former are in the 5 — 20% range, 
while for the latter, in the 20 — 80% range. As a re¬ 
sult, in one 0.25 s time bin, the flux typically changes by 
10 % after the true touchdown point, leading to a drop 
in the color correction factor from ~ 1.9 to 1.6, and a 
corresponding change in the blackbody normalization by 
50%. During the same time bin but prior to the touch¬ 
down point, the (i?app/D)^ changes also by ~ 50%, in 
the opposite direction. Therefore, it becomes impossible 
to detect the expected rapid evolution of the blackbody 
normalization as predicted by theoretical models. 


Looking for the evolution of the blackbody normaliza¬ 
tion as a function of flux fails not only on a burst-by- 
burst basis but also when it is applied to all the bursts 
from a given source. Fundamentally, the color correction 
factor for a given source depends on the local flux in the 
atmosphere, as measured by the effective temperature, 
and not on the total observed flux. The latter is affected 
by a number of factors, including the time bin averag¬ 
ing discussed above, as well as effects related to, e.g., 
partial obscuration of the neutron star surface or any 
reflection off of the accretion disk which are thought to 
be responsible for the outliers observed in some sources. 
This is of particular importance for the few sources that 
show a substantial number of outliers. The same effects 
also introduce a scatter to the blackbody normalization. 
Therefore, when the quantity is plotted as a func¬ 

tion of the measured flux, as shown in Figure |4] for the 
case of 4U 1608—52, the theoretically expected depen¬ 
dence is masked by the scatter. If we plot, instead, the 
quantity against the observed color temperature, 

the scatter in the abcissa is reduced to a point where the 
theoretically expected correlation starts becoming evi¬ 
dent. 

In order to demonstrate the degree of scatter that even 
a small spread in the emitting radius can introduce to the 
expected correlations, we performed the following Monte 
Carlo simulation. We simulated mo ck observations of 
burst spectra using Model #17 of ISuleimanov et al.l 
(j 2012 f l for the dependence of the color correction fac¬ 
tor on the effective temperature, which corresponds to 
a helium atmosphere with surface gravity log 5 = 14.3. 
We assigned a spread in the radius of the emitting region 
with a Gaussian likelihood. We used these mock obser¬ 
vations to infer the values of the blackbody normaliza¬ 
tion, color temperature, and flux that would have been 
measured under these conditions. Figure [5] shows the re¬ 
sulting relations between blackbody normalization, color 
temperature, and flux, for two cases: the left columns 
with a 5% spread and the right columns with a 15% 
spread in the emitting radius. Even for the small spread, 
the blackbody normalization appears a lot more tightly 
correlated with color temperature than with flux. Eor 
the larger spread, the latter correlation disappears and 
would have given the impression that the observations 
are inconsistent with theoretical expectations. Note that 
in performing these simulations, we have not included 
the significant smearing that occurs due to the size of 
the time bins that we discussed above or the statistical 
uncertainties that lead to correlated values between the 
blackbody normalization and the blackbody temperature 
(see, e.g.. Figs. 3 — 8 of lOzel et al.ll201^ . 

The above discussion demonstrates tha t the bursts 
selected in the data-driven approach of iGiiver et al.l 
(|2012airb| ) show a spectral evolution that do not contra¬ 
dict theoretical expectations and can, therefore, lead to 
an unbiased sample for the measurement of neutron star 
radii (as reported bv lOzel et al.ll2015h . Moreover, as we 
will show in the next section, devising selection crite¬ 
ria based on this theoretical expectation, without taking 
into account the inherent limitations of the data, leads to 
selecting precisely those bursts that do not reach the Ed¬ 
dington limit, therefore, violating the basic assumption 
of the procedure. 
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Fig. 5. — Monte Carlo simulations of the correlations between the blackbody normalization, color temperature, and flux based on a 
theoretical model of a bursting neutron star atmosphere. The theoretical model dependences, from which the mock observations were 
drawn, are shown by the red curves. In the left panels, we assigned a 5% Gaussian spread in the radius of the emitting region, whereas in 
the right panels, we assigned a 15% spread. Even for the small spread, the blackbody normalization appears a lot more tightly correlated 
with color temperature than with flux. For the larger spread, the latter correlation disappears and would have given the f^alse impression 
that the observations are inconsistent with theoretical expectations. 


3. PROBLEMS WITH THE RADIUS MEASUREMENTS IN 
THE THEORETICALLY SELECTED BURSTS 

iPoutanen et alJ (|2014D applied the criteria based on the 
theoretical evolution of the blackbody normalization to 
select bursts for the measurement of the neutron star ra¬ 
dius in 4U 1608—52. Specifically, they chose bursts where 
the blackbody normalization at half the Eddington flux 
is at least twice as large as the blackbody normalization 
at what they marked as the touchdown point in their 
sample of Eddington limited bursts. 

The identification of Eddington limited bursts relies 
exclusively on finding evidence for photospheric radius 
expansion, where there is a substantial increase in the in¬ 
ferred apparent radius of the photosphere above the value 
measured in the cooling tail, coincident with a decrease 
in the color temperature. Figure [S] shows the evolution of 
the spectral parameters in two bursts from 4U 1608—52. 
The burst on the left is rep r esent ative of all the bursts 
selected by IPoutanen et alJ (|2014ll as PRE bursts that 
follow the theoretical expectation. It is evident that even 
though the blackbody normalization and the color tem¬ 
perature oscillate in the beginning of the burst, the pho¬ 
tosphere never expanded, i.e., the blackbody normaliza¬ 
tion never exceeded the asymptotic value at the cooling 
tail. In contrast, the burst on the right is a bona fide 
PRE burst and shows the increase in the blackbody nor¬ 
malization during the radius expansion episode by more 
than a factor of thirty. Furthermore, the flux at the point 
that is identified as the touchdown on the left is substan¬ 
tially (> 30%) smaller than the Eddington flux measured 
at touchdown for the burst on the right. Even though 


the method followed bv IPoutanen et all (l2014ll relies on 
measuring the evolution of the color correction factor at 
flux levels within 10% of the Eddington limit, the selec¬ 
tion procedure failed to identify bursts that reach such 
flux levels. 

It is also worth emph asizing that, even though 
IPoutanen et al.l (1201411 and iKaiava et all (120141) charac¬ 
terize the selected bursts as long bursts occurring in the 
hard state, implying an underlying physical reason for 
this selection, their final sample represents neither the 
ones that occur in the hard state nor the longest ob¬ 
served bursts. This can be easily seen in Figure [7] that 
compares the distribution oi (i) the burst locations in the 
color-color diagrams and (ii) the time from burst start 
to touchdown between the selected bursts and the entire 
sampl e for 4U 1608—52. Figures 3 and 4 of iKaiava et al.l 
(|2014ll corroborate these observations and challenge the 
interpretation that the theoretically motivated burst se¬ 
lection has an underlying physical mechanism (e.g., one 
that is related to the different states of the accretion 
flow). 

4. CONCLUSIONS 

In this paper, we explored two approaches that have 
been proposed to select thermonuclear burst data for the 
spectroscopic measurements of n eutron star radii. T he 
first is the data-driven approach of lGfiver et al.l (j2012allbl l 
that used a Bayesian Gaussian mixture algorithm to 
identify data outliers and measure the degree of sys¬ 
tematic uncertainty in the measurement s. The second is 
the theoretically motivated approach of IPoutanen et akl 
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Fig. 6.— The evolution of the spectral parameters of two bursts from 4U 1608—52. The left panel shows one selected bv IPoutanen et alJ 
(|2Q14I') as a PRE burst, whereas the right panel shows a bona fide PRE burst according to lGiiver et al.l (|2012al) . The blackbody normalization 
in the first case never exceeds the asymptotic value at the cooling tail and its touchdown flux is substantially be low the touchdown flux of 
the second burst. This strongly argues against the classification of the bursts selected bv IPoutanen et al.l II2014I) as PRE bursts. 


(|2014D that required selected bursts to follow trends ex¬ 
pected &ran_the_buxsti^ neutro n star atmosphere mod¬ 
els of ISuleimanov et al.f (l2012fl . The two approaches 
led to non-over l appin g data selections. F urthe rmore, 
iPoutanen et al.l (|2014f l and iKaiava et al.l (|2014f l have 
used the theoretically motivated model to argue that the 
data-driven selection is not reliable because it is in con¬ 
flict with the theoretical expectations. 

We showed that the theoretically expected trends be¬ 
tween different spectroscopic quantities, such as the 
blackbody normalization, the color temperature, and the 
blackbody flux, get substantially smeared in the data be¬ 
cause of three effects. First, the inherent limitations of 
the RXTE data do not allow resolving the rapid spec¬ 
tral evolution that occurs at the end of the photospheric 
radius expansion episodes. Second, even a mild scat¬ 
ter in the emitting area, e.g., due to uneven burning, 
an evolving photosphere, partial obscuration of the sur¬ 
face, and/or reflection off of the accretion disk, lead to a 
much larger scatter in the measured spectroscopic quan¬ 
tities that mask the theoretically expected trends. Fi¬ 
nally, counting statistics lead to correlated measurement 
uncertainties between the spectroscopic quantities, fur¬ 
ther smearing any trends. We quantitatively showed that 
these effects are at a level to prevent the detection of 
theoretical trends in the majority of the burst data and 
preclude using these criteria to argue against the data- 
driven selection methods. 

Finally, we studied the sub-sample of bursts selected 


bv IPoutanen et all (j2ni4[ l for 4U 1608—52 using the the¬ 
oretically motivated criteria without taking into account 
the data limitations discussed above. We showed that, 
contrary to the implicit assumption in the method, the 
bursts selected are not PRE bursts. Indeed, at no point 
during any of these bursts does the radius of the pho¬ 
tosphere exceed the asymptotic radius measured during 
the cooling tails of the same bursts. Furthermore, the in¬ 
ferred touchdown fluxes of these bursts is 30% below the 
touchdown fluxes of the securely identified PRE bursts 
that reach the Eddington limit. 

The quantitative analyses presented in this paper leads 
us to the conclusion that this theoretically motivated 
selection procedure is neither practical, given the limi¬ 
tations of the data, nor unbiased, given that it selects 
bursts that are inconsistent with its own assumptions. 
However, in the future, with the use of an X-ray detec¬ 
tor with a larger collecting area and significantly fewer 
limitations when observing sources with high count rates, 
this could lead to useful radius measurements. 
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Location in Color-Color Diagram S2 Characteristic Timescale (s) 

Fig. 7.— The distributions (Left) of hardness in the persistent emission (as measured by the location S-z along the color-color diagram) 
and (Right) of the time between burst start and touchdown, for (blue) all the bursts of 4U 1608—52 analyzed by Poutanen et al. (2014) 
and (red) for those selected by them as being consistent with the theoretical predictions. The selected bursts are neither the longer bursts 
nor the ones that occur in the hardest state, contrary to earlier claims. 
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